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We investigate the quantum roll for a particle in a d-dimensional “Mexican hat” potential 
in quantum mechanics, comparing numerical simulations in d-dimensions with the results 
of a large-d expansion, up to order 1/d, of the coupled closed time path (CTP) Green’s 
function equations, as well as to a post-Gaussian variational approximation in d-dimensions. 
The quantum roll problem for a set of N coupled oscillators is equivalent to a (d = N)- 
dimensional spherically symmetric quantum mechanics problem. For this problem the large- 
N expansion is equivalent to an expansion in 1/d where d is the number of dimensions. We 
use the Schwinger-Mahanthappa-Keldysh GTP formalism to determine the causal update 
equations to order 1/d. We also study the quantum fluctuations (r^) as a function of time 
and find that the 1 /d corrections improve the agreement with numerical simulations at short 
times (over one or two oscillations) but beyond two oscillations, the approximation fails 
to correspond to a positive probability function. Using numerical methods, we also study 
how the long time behavior of the motion changes from its asymptotic (d ^ oo) harmonic 
behavior as we reduce d. 
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I. INTRODUCTION 

The quantum roll problem is the problem of studying 
the quantum behavior of a particle starting in an unsta¬ 
ble equilibrium at the top of a potential hill and “rolling” 
down with the constraint (x) = 0. The roll problem re¬ 
cently became relevant to cosmology with the advent of 
the new inflation model [Q. In that model of the early 
universe one is interested in determining the slow rollover 
of the scalar inflaton field in a time evolving semiclassi- 
cal gravitational field. The quantum mechanics of the 
slow rollover has been previously studied mostly in per¬ 
turbation theory or in a mean field approximation (see 
however [^). For the one dimensional roll problem, it was 
found that the mean field approximation broke down be¬ 
fore the particle reached the bottom of the well. Thus it 
is important to find approximations which are valid for 
longer periods of time. One such approximation is to in¬ 
clude the 1/A^ corrections to the mean field result. Our 
interest in the quantum mechanical model was due to the 
fact that for this particular problem one is able to per¬ 
form numerical simulations at arbitrary N. This is due 
to the radial symmetry of the roll problem. If we were 


‘electronic mail:bogdan.mihaila@unh.edu 
^electronic maihjohn.dawson@unh.edu 
■^electronic mail:fcooper@lanl.gov 
^electronic mail:me_brew@c3serve.c3.lanl.gov 
“electronic mail:habib@lanl.gov 


interested instead in studying tunnelling numerically in 
an N dimensional anharmonic oscillator, one would be 
restricted to < 5, using even the most powerful com¬ 
puters available 1^. 

This paper is constructed as follows: In section || 
we state the d-dimensional quantum roll problem and 
discuss how to choose initial conditions which a llow a 
smooth transition to the N oo limit. In section III we 
derive a Post-Gaussian variational approximation valid 
at arbitrary N. In section IV we review the update equa¬ 
tions for the Green’s functions valid up to order 1 /N. In 
section ^ we compare both approximations to exact nu¬ 
merical simulations. In section VI we study the long time 
behavior of the two-point function as a function of N. 
We reserve the appendices for discussing our numerical 
methods. 


II. QUANTUM ROLL 

The quantum roll problem for an upside down Har¬ 
monic oscillator with Hamiltonian H = (p^ — a;^)/2 is 
the starting point for our discussion. The solution for the 
spreading of the wave packet with the constraint {x) = 0 
is exactly known Q. For example, if at t = 0 the wave 
function is Gaussian, 'if{x) = exp{—a:^/2}, then for t > 0, 

{x'^(t)) = 1/2-1- sinh^(t) . 

This exponential growth of course does not continue in a 
bounded potential. In this work we will study the quan¬ 
tum roll for a set of N coupled oscillators in a “Mexican 
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hat” potential where the motion is bounded, as a function 
of the number of oscillators. We take the Hamiltonian to 
be of the form, 

H = -^y^ + V{r), ( 2 . 1 ) 


where 


= "ils + t'M. P-9) 

with an effective one dimensional potential U(r) given by 


where 


N 


v^ = E 


2=1 


dx? ’ 


N 




r(P = (P - rl) 


2=1 
2\2 


with the Schrhdinger equation given by: 


H = i 




dt 


( 2 . 2 ) 


(2.3) 


It is convenient to study this problem in a multidimen¬ 
sional coordinate system with d = N, using the radial 
coordinate r and a set of d — 1 angular coordinates. We 
write the Laplacian as 0,^: 


E , (^- 1 ) 9 _ 

dr '2 r dr 7.2 ’ 


(2.4) 


where is the generalized orbital angular momen¬ 

tum operator defined by Louck Q . We separate the wave 
function as as radial function times a hyperspherical har¬ 
monic: 


^(Xi,t) = 4'djA](»',i)^d.[A](f^) , 
then the equation for 4'^ [;^](r, t) becomes: 

9^'d.[A](?', t) 


H{r,l)'^d,[\]{r,t) = i 


where 




1 f d—ld 


2 \dr^ 


dr 


dt 


l{l + d-2) 


(2.5) 


( 2 . 6 ) 


V{r) . 


Here A denotes the set of (d — l)quantum numbers, in¬ 
cluding the orbital quantum number I needed to specify 
completely the hyperspherical harmonic (see for exam¬ 
ple g). 

F or t he quantum roll problem the potential is given 
by (|2.2| ). If we start at the top of the hill at r = 0 
with a radially symmetric initial state centered at r = 0 
then there is no angular momentum and we only have to 
consider the case I = 0. The first order derivative term in 
the Hamiltonian can be eliminated by the substitution, 


U{r) 


(d-l)(d-3) 

8 r 2 



riy . ( 2 . 10 ) 


Using this Hamiltonian we can update $ using simplectic 
methods, or by solving numerically for the eigenfunctions 
and eigenvalues and using the expansion: 


^(r,t) = r-^d-i)/2 ^ c„e-*^"‘$„(r), (2.11) 

n 


with Cn being determined from the initial conditions on 
4^ using the orthonormality of the 4)„(r). 


A. Initial conditions 


For most of this paper, we will use for simplicity, Gaus¬ 
sian initial conditions, since they allow for a simple deter¬ 
mination of the Cn- However, in this paper we also want 
to consider a generalization of the Gaussian variational 
method discussed in g, which was found quite robust 
in describing the time evolution of pulses in classical dy¬ 
namical systems. This generalized wave function has the 
ability to approximate in shape an arbitrary spherically 
symmetric pulse that is monotonically decreasing around 
r = 0. This is accomplished by adding one more varia¬ 
tional parameter a which changes the shape of the pulse 
from flat to peaked, including the Gaussian as the spe¬ 
cial case a = 1. Pulses of this form stay coherent for 
a long time in many nonlinear equations and have been 
used previously in variational calculations of soliton mo¬ 
tion to study soliton blowup. They have the advantage of 
allowing analytical expressions for the expectation values 
and thus are ideal for the variational calculation we will 
consider below. 

Thus we will consider at < = 0, normalized wave func¬ 
tions of the form 


^h(uO) = 


2a 


- 11/2 


ndT[d/i2a)] 
X exp 


2a G 


[2a] 


-d/{4a) 


4aG| 


[2a] 


( 2 . 12 ) 


where fid is the angular volume, = 2 r(d/2), 

such that we satisfy the normalization condition 


$(r,f) = r(‘'-i)/2 4r(r,t) , (2.7) 

in which case, the time dependent Schrodinger equation 
for $ becomes: 

H\r,l) ^r,t) = z , (2.8) 


1 = 


poo 

/ |^[a](uO)P dr. (2.13) 

^0 


The energy of this state has the form: 

Eia) = a(a)Gf4“ + e(a) gJ/^^ + /(a)G;/“ . (2.14) 
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where 


a 


2a 


— 1/a. 


a{a) =-{d + 2a-2) ( — ] R{2a-2,a), (2.15) 


s(a) = 


8d \ d ) 


2/a 


R{A,a) , 


l/a 


and we have defined i?(/3, a) by: 

/? + d 


Rif3,a) = r 


2 q! 


, , d 


(2.16) 

(2.17) 

(2.18) 


In particular, for a = 1, the above initial conditions be¬ 
come 


^[ 1 ](^. 0 ) = 


1/2 


2G 


[ 2 ] 


_r!dr(d/2)J 
We notice that Gild = G where 
G = {x\) 


-d/4 — 


4G| 


[ 2 ] 


for each i. Evaluating the energy (2.14) for the Gaussian 
initial state {a = 1), we obtain: 


E=^ [d{d + 2)G^ -2rldG + r^l + 

od 


8 G' 


(2.19) 


B. Asymptotic form of the wavefunction for large d 

The eigenfunctions in the large-d limit are Gaussians 
times polynomials, where the Gaussians are centered 
about the minimum of the effective potential of the one- 
dimensonal radial problem. In order to get a uniform 
overlap at arbitrary d it is important to choose initial 
conditions so that the expansion in terms of eigenfunc¬ 
tions is similar at all d. This requires us to run the cou¬ 
plings as a function of d so that the overlap is constant 
up to terms of order 1/d^. This can be done in several 
ways that differ by terms of order 1/d^. The method 
presented below leads to uniform results even at d = 1 
as we change the parameters with d. 

In order to do this we need to examine the asymptotic 
form of the wave function in the large d limit, and choose 
parameters so that in this limit, we do not introduce 
undesirable numerical errors in the initial decomposition 
of the wavefunction into eigenmodes. At t = 0, the initial 
wave function 4'o(r’) is a normalized gaussian centered 
about the origin: 


\I'o(»’) = V exp 



with normalization Af, such that Eq. ( 2.13| ) is satisfied. 
However the rescaled wave function 4>o(r') at t = 0 is 
given by: 


$o(r) = exp 


= Vexpi-- 



Thus the rescaled wave function, for large d, can be ap¬ 
proximated by a gaussian centered about r = Too- That 
is, for large d. 


$o(0 « Af exp + 0(l/d) )■ , (2.20) 


which defines foo and Goo- For all d, however, we have: 

nOC 

f = (r)^ ^ 

0 


= 


F((d+l)/2) 


r(d/2) 

pOO 

J 0 

= VdG. 


Thus, we define: 
Goo 


( 2 . 21 ) 


( 2 . 22 ) 


- {ry = G<d-2 


r((d + i)/ 2) 
r(d/2) 


1 2 


or 


G = 


Go 


2d-4 


r((d + i)/2) 


(2.23) 


r(d/2) 

In the limit when d goes to infinity, we have 
G(d) 


Go 


f(d) 


Too = V{d- l)Go 


Now the solution 4)(r, t) satisfies Eq. ( 2.8) with the po¬ 
tential function U (r) given by Eq. ( ^.IC ). In the spectral 
method, we expand solutions of Schrodinger’s equation 
in eigenvectors of this equation. Eor large d, these eigen¬ 
vectors are centered about the minimum of the potential 
I7(r). U{r) has the expansion. 


U{r) = U{r) + -w? {r — ry + ■ ■ ■ , 
where f is given by the solution of the equation, 


(d-l)(d- 3 ) 

4 


9 r-2 


(2.24) 


(2.25) 


and ra^ by: 
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w? 


3{d l)((i 3) g , n 2 \ 

4f4 ° ’ 

^ (3f2 - 2rg) . 


(2.26) 


We want to make sure that the difference Sr between r 
and r, 


Sr = r — r , 


remain a constant for all d. This will insure that the 
overlap integrals for initial values of the coefficients in 
a secular expansion will remain approximately the same 
for all values of d, and thus numerical errors associated 
with the initial conditions will not effect our results. 

We also define to be the second derivative of V (r) 
evaluated at r = f, 


2 _ d‘^V{r) 


then, 


dr'^ 
becomes: 




\d)=m^- 


_2 3(d-l)(d-3) 


4r 




(2.27) 


(2.28) 


Solving ( ^.27l ) for ro, substituting into ( 2.25| ) and solving 
for g gives: 


aid) 



id-l)id-i) \ 

^4 J 


(2.29) 


Holding ifi fixed means that the frequency of the oscilla¬ 
tion remains the same for all d. Holding Sr fixed, means 
that the overlap between the initial wave function and ini¬ 
tial starting values for the solution remains consta nt fo r 
all d. T hus fo r fixed value s of Sr, Goo and Eqs. ( ^.21 ), 
(2.23), ( 2.28| ), and ( ^.29 ), determine values for G, f, f, 
and g for all values of d. Thus in order to keep the 
same accuracy in our solutions, we need to run the cou¬ 
pling constant g{d) and the initial width G{d) as defined 
above. 

We define poid) and p{d) by the ratios. 


Paid) = rl/{Gd) , p{d) = f'^/{Gd) , 


where Tq is calculated from: 


2 -2 
Tn = r 




1 - 


-2-4 

mr 


(d- l)(d-3) 


(2.30) 


Note that in the limit when d goes to infinity, the various 
parameters discussed in this section have the limit 


aid) - 

aioo) = 

1 


(2.31) 

mfi{d) 

TO^(oo) = 

- 

3 

4G2 

(2.32) 

paid) - 

> po(oo) = 

I -f 

I 

(2.33) 

2 (l-m 2 G^) 

aid) 

^ p{oo) = 

I 


(2.34) 


We show in Fig.|| a plot of po(d), /5(d), wfi{d), and g{d) 
for the cases Goo = 1 , dr = 2 , and md = 2 and m? = 
21. The corresponding initial values are given in Table || 
and ||. These plots show how these parameters flow as 
a function of d. This insures that the initial Gaussian 
wave function would have finite overlap with the infinite 
d eigenfunctions. Other choices which differ by terms of 
order 1 /d^ are also possible. 

If we are using the spectral method for solving the 
Schrodinger equation, we need to make sure that the en¬ 
ergy of our initial wave packet is not too large compared 
with the ground state, so that the solution will be valid 
for long times. In order to estimate the ground state 
energy it is sufficient to make a harmonic approximation 
around the minimum of th e effec tive one dimensional po¬ 
tential U{r), given by Eq. ( ^.24 ). The ground state wave 
function is 


$(r) = N exp --^(r — r)^ 


(2.35) 


where is a normalization constant. Correspondingly, 
the ground state energy is approximately 


Eq = Uir) -I- m/2. 


(2.36) 


We used this result to get a feeling for how to choose our 
initial conditions. Afterward we determined our parame¬ 
ters from the asymptotic relations, keeping this result in 
mind. 


III. TIME-DEPENDENT POST-GAUSSIAN 
VARIATIONAL METHOD 

Because of the radial symmetry of the quantum roll 
problem, it is possible to make a simple generalization 
of the Gaussian approximation which allows us to track 
wave packets which are have arbitrarily high correlation 
functions and thus might be more coherent in their be¬ 
havior in a non-linear potential. These wave packets were 
used previously to study soliton behavior and blowup 
in the classical nonlinear Schrodinger equation in arbi¬ 
trary d. These wave functions have the property that at 
d = 1 they give a much better ground state wave func¬ 
tion in that the wave function becomes quite flat and can 
span both wells. At large d this wave function goes over 
to the gaussian limit since that is the correct ground state 
at large d. Since the analytic large-d expansion is an ex¬ 
pansion about a Gaussian, we expect the post-Gaussian 
variational wave to be compatible with that expansion 
only at large d. At small d, if we use the parameters ob¬ 
tained from the ground state energy, our initial state will 
be quite different from a Gaussian and not track an initial 
Gaussian very well in detail. We can, however start with 
non-Gaussian initial conditions and make direct compar¬ 
isons numerically. 

The time dependent Schrodinger equation can be ob¬ 
tained by varying the Dirac action 


4 
































r = J dt{'i/{t)\id/dt-H\'i/{t)) . (3.1) 

The Schrodinger equation results by asserting that the 
action T be stationary against arbitrary variations of the 
wave function 'I' (r, t). We approximate the true wave 
function (r, t) using a generalized Gaussian trial func¬ 
tion of [|| 


'^v{r,t) = 


2a 


ndT[d/{2a)] 

,2a 


X exp 




— iA 


4aGp 


Correspondingly, the expectation value of is 

„ /2aGr2ai\^/^^“^ 

{r^)=[ --—I) R{f},a), 


(3.2) 


(3.3) 


where R{f3,a) is given by Eq. (p.l^). In particular for 

/3 = 2a, 

G[2a] = (r^“) . 

Note that in ref. j|], we considered only the case when 
a = I where 

G[2] = (r2) = dG. 


Using Eq. (^), we can calculate the expectation value 
of the Hamiltonian, {H) = (iJfree) + {Hint), with 


(-^free) — ^ ^ 

u 


i-l.T,* 




d-1 d 


dj.2 


r dr 


4'^(r, t) 


= a {d + 2a — 2) 

r, . . d 

-2a^ 


^4aG 
= a {d 2a — 2) 


4aG 
— iA ] (r‘ 


- lA ' 


4aG 


-A" 


m 


T 


(2a-2,a) 


(3.4) 


and 


(ffin.) = njf fUr.t) V(r) <t„(r,t)r‘-'dr 
We also find the expectation value of i d/dt as 


{i d/dt) ^ 


^11 ^ . ^7. 


dt dt 


— ~ A , 


(3.6) 


so the effective action can be written as: 

r = y'dt {-G[ 2 „]A - a(a)G[-4“ 

- - e(a) G^/“ - f{a) G^^ } , (3.7) 

where the coefS c ients a(a), e{a) and /(a) are given by 
Eqs. (2.15, 2.16, 2.17|). In addition we have introduced 


b{a) = a(a) 


4a 


(3.8) 


The variational equations of motion are obtained from 


M 


= G[2„] - 26(a)G^4/“A = 0, 


(3.9) 


sr 


SG 


= -A- 


[2a] 


a(a) i-i/c 
’^I2a] 


a 


-6(a) 2-- G 


a 




b2a] 


2e(a)^_i+2/c 

^[2a] 


/(«J ^-1 + 1/a 


+^/“ = 0. (3.10) 


a ‘ J a 

Note that Eq. (^) can be solved for A. We get 


A = 


1 ^-2-l-l/a 


(3.11) 


Since the effective action is stationary against variation 
of the tria l wave function 'I'^, (t), the equations of motion 
( |3.9i tl.lO ) guarantee energy conservation: 


Eld = a{a)G-^i{^ + 6(a) Gf -4 


'(2-l/a^2 


+ e(a)G^2,^ + fi^)G[2a] ■ (3-12) 

The conservation provides a first integral of the motion 
and gives us a non-linear differential equation to solve in 
terms of G[2a] and the initial energy. 

Stationary solutions are obtained by requiring G[2a] = 
A = 0. Then, from Eqs. (3.9, ^.10 ) we obtain A = 0 and 
an equation for Gpa] 


yl/a 


2e(a)Gf/“ + /(«) - a(«) = 0. 

The energy of this stationary solution is: 


(3.13) 


E{a) 


— a(a)Gj2^^“ + e(a)Gj2''“] + /(a)Gp^ . 


2/a 


.,1/a 


By minimizing this energy with respect to G[2a] and a we 
get the best estimate of the ground state wave function 
for this class of trial wave functions. We must remember 
that we expect the actual ground state wave function to 
have its support at f and not around zero. This gen¬ 
eralized approximation has the feature that as d 00 , 
ac 1 . 

First let us look at the ground state wave function 
at d = 1, (?oo = 1) which, from Table 0, corresponds 
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to g = 0.232 and tq = 2.9359. We compare the vari¬ 
ational wave function and ground state energy with the 
best Gaussian result. Minimizing the energy with respect 
to the parameters G[ 2 a] ^.nd a, we find that the optimal 
parameters are: 

ac = 5.25774 ; G[2a,] = 152039., (3.14) 

leading to a normalized wave function given by: 

V'(r) = 0.134818 exp {-6.25428 x j-io.5i5483 | ^ 

corresponding to an energy of E{ac) = —1.05847. In 
comparison, the Gaussian approximation (a = 1) leads 
to the results: 

a = 1.0; G[2] = 2.95542 , (3.15) 


Gab (^> ^') = I ^ + X(i) I <^c (t, t') Sab 

= GoHt,t')dab, (4.2) 

= (4.3) 

N 

Gab{t,t')Gba{t',t) . (4.4) 

a,6—1 

Here 6c{t,t') is the closed time path delta function. The 
auxilliary variable x(0 obeys the constraint equation 
given by 

A ^ 1 

X(t) =-p — ^-^„„(t,t) , (4.5) 

a—1 


with a wave function given by: 

V'(r) = 0.232058 exp {-0.16918 r^} 

corresponding to an energy of E{a = 1) = —0.67531. A 
comparison of these two wave functions is found in Fig. |^. 
Thus the more general wave function has support over 
both wells and does considerably better for the energy. 

Because of this wider support, the post-Gaussian ap¬ 
proximation does a better job in reproducing the am¬ 
plitude of the excursions in the time dependent prob¬ 
lem when compared to the Gaussian approximation even 
though it does not represent the early time evolution. 
At larger d > 10 when a becomes closer to unity then 
the entire evolution is modelled better by this type of 
trial wave function. This is seen in Fig. ||. We can also 
compare the generalized Gaussian trial wave function re¬ 
sult with with a numerical simulation starting from the 
non-Gaussian initial conditions using the split operator 
technique described in . The results at two different d’s 
are shown in Figs. Again we see that we need d > 10 
for this approximation to give reasonable results. 


IV. LARGE-A EXPANSION 

The large-A (d) expansion for the d-dimensional an- 
harmonic oscillator is the zero space dimensional limit of 
the (j)^ field theory formalism discussed in ref. §. The 
effective action to order 1/d is given by [here we use the 
field theory notation of ref. g: d ^ N, g ^ X, and 
^2 ^ -grl/{2d)], 


and where the full Xi propagator and self energy 

to order 1/iV are given by 

Qab{t,t') = Gah{t,t') (4.6) 

~ Yh [ / <^^2 Gac(tj G) Scd(tl) ^ 2 ) G£;b(t2) ^0 

Je 

Scd(t, t') = i Gcd{t, t') D{t, t'). 

These equations are also derived in (2.18-2.22) of ref. [||. 
In order to solve for we first write 


= -6cit,t') + j AD{t,t') . (4.7) 

Then we find that AD{t, t') satisfies the integral equa¬ 
tion. 


AD{t,t') = An(f,t')-^ciGn(f,G)AG(t",d), 

(4.8) 

in agreement with (2.13-2.16) of Ref. |^. We are now in 
a position to solve these coupled equations for the motion 
of GabjEt') and x{t) for given initial conditions. We solve 
Eq. (|/^) simultaniously with (4.6) and (4.8), using the 
Ghebyshev expansion technique of appendices A and B 
of ref. [0. 


V. COMPARISON OF APPROXIMATIONS WITH 
EXACT NUMERICAL SIMULATIONS. 


r[x,a;,] = 

/ “ XxUt)] + ^ln[Go^(<,t)] 

c ,• 


+ 


Nx^it) Ng'^xit) , * 


2A 


A 




(4.1) 


where G^^{t,t') and D are the inverse propaga¬ 

tors for Xi and x, given ( for (xi) = 0) by 


In this section we will compare the exact solution (de¬ 
termined numerically) with both the variational method 
and the large-A expansion. In order to have a smooth 
d —> 00 transition, we will allow our parameters to run 
as a function of d as described above. We use the set of 
parameters listed in tables | and ||, which correspond to 
the choice of Goo = 1, Sr = 2, and the two cases fri^ = 2, 
(7oo = 1 and fn^ = 21; g^c = 20, respectively. First let us 
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discuss Gaussian initial conditions. In Figs. we present 
the results for Lr^ jd) which correspond to goo = 1. We 
find that at short times (less than 2 oscillations) the cal¬ 
culation which includes 1 jN corrections tracks the exact 
result the closest, however when it starts deviating from 
the exact result it can lead to a negative expectation 
value. This is evidence that the next to leading order 
large-N approximation does not correspond to a positive 
definite probability function. At leading order, the large- 
N approximation is equivalent to a Gaussian wave func¬ 
tion (or density matrix) and all the expectation values 
of even moments of operators are positive. This problem 
of not having a positive definite probability associated 
with this type of approximation scheme also arises if we 
consider approximations which are exact truncations of 
the Green’s function hierarchy |10| at the level of the 
connected four point function. 

The leading order in large-N calculation is closer in 
amplitude to the exact result than the Hartree varia¬ 
tional approximation. These latter two approximations 
do not suffer from the illness of the next to leading large- 
N approximation. We see that all the approximations 
approach the exact one as d ^ oo as they must. In 
Figs. ^ we see the same type of behavior at a larger value 
of goo = 20. The main difference here is that the time 
scale for a single oscillation has been reduced, but the 
behavior in terms of the oscillation time scale is similar 
to the previous case. 

We have also compared the post-Gaussian variational 
approximation for Gaussian initial data.This comparison 
is shown in Figs. ||. Here, the initial states are clearly dif¬ 
ferent so we have normalized the two wave functions so as 
to have the same value for (r^). What is interesting here 
is that the post-Gaussian approximation gives a better 
approximation for the amplitude of the oscillation, which 
’’cures” the previously found failure of the Gaussian ap¬ 
proximation to reach the bottom of the Mexican Hat po¬ 
tential in one dimension. To truly test the post-Gaussian 
approximation we used non-Gaussian initial conditions 
pertaining to the value of a = Uc found by minimizing 
the ground state energy. To obtain our numerical results 
here, it was not possible to use the eigenfunction method 
because of the number of eigenfunctions required to ap¬ 
proximate this type of wave function. Instead we used 
the split operator technique of ||^. As shown in Figs. Q we 
see that the post-Gaussian approximation does quite well 
with the amplitude of the oscillation, but only at large 
d > 10 does it start getting the frequency of oscillations 
correct. 


VI. LONG-TIME BEHAVIOR OF (r=) AS A 
FUNCTION OF d. 

At large d the effective one dimensional radial poten¬ 
tial is very deep and becomes harmonic for oscillations 
around the minimum. Thus as d —> oo an initial Gaus¬ 


sian wave packet should stay gaussian, because the gaus- 
sians are coherent states of the harmonic oscillator. How¬ 
ever, the solutions for finite dimension d can not be ex¬ 
pected to maintain the coherence of the original Gaussian 
state over long time intervals. It is of interest to observe 
the actual long-time behavior of the solutions that arise 
from Gaussian initial conditions when the potential is not 
purely quadratic. 

The method of eigenfunction expansion can be used 
to calculate the solution over large time intervals. The 
accuracy of the long-term behavior is limited primarily 
by the accuracy of the eigenvalues, so the tolerance of 
the eigenvalue calculation was reduced until numerical 
convergence of the solution was observed over the time 
intervals of interest. 

The few examples that we consider in this paper show 
a surprising variety of long-time behaviors. In Figs. 
we show the evolution of G 2 = (’’^/d) over a time interval 
sufficiently long to capture the characteristic behavior of 
that particular solution. 

In Fig. ^ we notice that as we increase d we gradually 
approach the limiting form of oscillation in an harmonic 
well. In the case d = 1, = 2, (Fig. ^ G 2 is domi¬ 

nated by the contribution of the lowest frequency and so 
appears as a perturbation of a simple oscillation with fre¬ 
quency equal to the energy difference between the ground 
and first excited states. The solution itself contains one 
dominant peak with a few smaller peaks. As we increase 
d the potential gets more harmonic and simple harmonic 
motion of G 2 is observed. At larger g^o, as seen in Fig. ||, 
more complicated behavior can be seen. 

In the case d = 1, = 21, (Fig. ^) there is a sig¬ 

nificant modulation of the amplitude, so that at regular 
intervals there is a much reduced oscillation of G 2 . The 
solution typically has several (2-3) peaks of similar size. 

In the other cases in Fig. |[ there is always significant 
amplitude modulation so that there may be brief or in 
some cases (d = 10 and 20, = 21), extended periods 

when the amplitude of the oscillation is an order of mag¬ 
nitude less than the initial oscillation. The oscillation 
always returns to a significant fraction (50-90%) of the 
initial oscillation, even after extended periods of reduced 
amplitude. The pattern of amplitude modulation can be 
quite complex. The solutions in general have a number 
of major peaks, from a minimum of 3 up to 12, indicat¬ 
ing the contributions from higher eigenstates. However, 
in some cases the coherence of the initial Gaussian state 
may be maintained for a significant length of time. For 
example, in the case d = 20, frd = 2, the coherent state is 
maintained approximately for 10 oscillations at the fun¬ 
damental frequency. For this large value of g we expect 
that one needs to get to very large values of d > 100 
before the asymptotic behavior seen in Fig. will be 
reached. 
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VII. CONCLUSIONS 


APPENDIX A: NUMERICAL APPROACH 


In this paper we have studied the behavior of G 2 = 
{r^/d) for Gaussian and non-Gaussian initial data for the 
quantum roll in d dimensions. We found that as d —*■ 00, 
G 2 became more harmonic in behavior as expected. We 
compared the exact numerical data with two variational 
approximations (Gaussian and Post Gaussian) for the 
wave function as well as the leading and next to lead¬ 
ing order approximations in the large d expansion for the 
Green’s functions. We found that all these approxima¬ 
tions converged to the exact result as d ^ 00 . At short 
times (less than 2 oscillations) the next to leading order 
large d approximation was the most accurate one. How¬ 
ever when this approximation started diverging from the 
exact result, it broke down in a serious fashion in that 
G 2 became negative. This unexpected result is related 
to the fact the large d expansion for the expectation val¬ 
ues may not necessarily correspond to a positive definite 
probability when truncated at any finite order in 1 /d. At 
lowest order, however the approximation is equivalent to 
a Gaussian density matrix and does not have this defect. 
Variational approximations to the wave function by their 
very nature never have this defect. Recently |]^ we have 
seen that in a related approximation, namely truncating 
the equal times connected Green’s function at the four 
point function level, G 2 is not positive definite. Thus we 
feel that truncating the heirarchy of Green’s functions 
in various ways might always have this problem, apart 
from the Gaussian case. This leads us to suspect that 
a higher order variational approach will be necessary if 
we want to include scattering corrections to mean field 
theory and also look at long time behavior. This will be 
the subject of a future paper. For astrophysical applica¬ 
tions of the quantum roll problem, one is only interested 
in getting the correct result to the bottom of the poten¬ 
tial, after which particle production leads to dissipative 
effects as one oscillates at the correct new minimum. For 
this purpose, keeping the next to leading order in large 
N definitely improves the mean field result. 
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1. Analytical Preliminaries 


The time-dependent Schrodinger equation 


19 ^ .d4> 

=‘a 


where U (r) is given by iH), 


U{r) = 


{d I)(d 3) 9/2 2\2 


8r2 




(Al) 


(A2) 


is solved by the method of eigenfunction expansion. Tak¬ 
ing a Fourier transform in time, we have the eigenvalue 
problem 


-^$" + C/$ = F;$ (A3) 

on the half-line r > 0. The boundary conditions are 
finiteness at the origin and the normalization condition, 

|$(r)|2dr=l. (A4) 

The numerical solution of this eigenvalue problem is de¬ 
scribed in the next section. The value of the wave func¬ 
tion at t = 0 is given by: 



4>{r, 0) = ipoir) = exp(—r^/4G) , (A5) 

where A is a normalization constant depending on G and 
d. So the initial-value problem is solved by the expansion 

00 

i) = Yl exp(-iF;jt)$j (r) (A6) 

j=o 

where the coefficients are determined by quadrature 

pOC 

Cj= Mr)<^j{r)dr. (A7) 

Jo 

To evaluate a particular moment of the solution without 
determining the full spatial solution, we take the moment 
of the expansion. Thus we may calculate the second mo¬ 
ment as follows 


= X! ‘^j,kCjCk cos{Ej - Ek)t , (A8) 

j,k=0 


where 


poo 

j,k= / r^$j(r)$fc(r)dr . 

Jo 


(A9) 




2. Numerical Solution of Eigenvalue Problem 

a. Discretization 


The ordinary differential equation (2) is discretized us¬ 
ing a compact finite difference scheme of second-order 
accuracy. This scheme allowed the easy use of a nonuni¬ 
form mesh near the origin, which was needed especially 
in the lower-dimensional cases {d < 10) to acheive the 
desired accuracy. The nonuniform grid {vn} was created 
by transforming the independent variable by the formula 




(AlO) 


where b = 2/{d — 3) and then using a uniform grid [h, 
2h, ..., L] in the variable s. The discretized equations 
are then 

ttn+l '^n _ '^n-t-1 T Wn 

^n+1 2 


and 


1 / Wn+l -Wn \ 

2 \ ^n-t-1 / 


^ Uir^+i) + Uir„) 


- E 


M(r„+i) +u{rn) 


= 0 . 


The boundary conditions at the singular points 0 and oo 
were imposed by specifying w = 0 at the leftmost grid 
point and rt = 0 at the rightmost grid point. The param¬ 
eters h and L were adjusted until sufficient accuracy was 
acheived, as described in the next section. 


b. Eigensolver 


In matrix notation we have Ax = E B x, where a; is a 
vector containing the u and w values at the grid points. 
This is a generalized eigenvalue problem and the matrix 
B is singular. We can formally convert this to a regular 
eigenvalue by inverting A, (although we will never actu¬ 
ally invert A in practice.) We then have A~^Bx = Xx, 
where A = E~^. This eigenvalue problem was solved for 
a specified number (16 usually) of the largest eigenvalues 
by the method of Arnold! factorization which is an 
iterative method and so is very fast provided matrix vec¬ 
tor multiplication can be carried out quickly. MATLAB 
routines where used for the Arnold! factorization as well 
as an incomplete L-U factorization of the matrix A that 
was used to quickly evaluate the matrix-vector products 
A~^Bx. 

Convergence of the solution was confirmed by testing 
the orthogonality of the eigenvectors and by computing 
the energy of the eigenfunctions from the formula 


E = 


rl- 

5$(r) 

Jo js 

dr 


+ V{r) |d)(r)p > dr , 


The size of the interval, as determined by L, and the 
mesh spacing, as determined by h, that were required to 
acheive an accuracy of 10“^ was highly dependent on the 
number of eigenvalues being determined and to a lesser 
extent on the dimension d. At most 4096 grid points 
were required in the worst case of 32 eigenfunctions with 
d = 20 — usually 1024 grid points were sufficient. 


3. Solving the Initial-Value Problem 


To solve the initial-value pro blem , the only remaining 
step is to evaluate the integrals (^). This is done by us¬ 
ing the trapezoidal rule. The accuracy of the quadrature 
was verified by reconstructing the initial condition and 
the mesh refin ed u ntil an accuracy of 10“^ was acheived. 
The integrals (^) were evaluated similarly. 
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0 1 2 3 4 5 6 t 

FIG. 2. Gaussian versus generalized Gaussian wave func¬ 
tions, for Qaa = 1 and d = 1 with ro = 2.9359. d=2Q 


FIG. 3. Gaussian initial conditions with goa = 1. Here we 
compare the exact results with both the post-Gaussian aprox- 
imation normalized to (r^) and the Gaussian approximation. 
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(a) d=l 



0 2 4 6 8 10 12 14 


(b) d=5 



(d) d=20 

FIG. 4. Post-Gaussian initial conditions with Qaa = 1. 
Here we compare the exact results with the Post Gaussian 
approximation with Oc fixed by minimizing the ground state 
energy. 
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(a) ci=l 
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(b) d=5 



0 5 10 15 20 

(d) d=20 

FIG. 5. Gaussian initial conditions with goo = 1. Here we 
compare exact result with the leading and next to leading 
order large-A^ approximation and to the Gaussian approxi¬ 
mation. 
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(d) d=20 


(d) d=20 


FIG. 6. Gaussian initial conditions with g^o = 20. Similar FIG. 7. Long-time behavior of the exact solution for Gaus- 

to Fig. 5 we notice that the next to leading order large-A'^ sian initial conditions with Qac, = 1, = 2, Goo = 1. 

approximation becomes unstable for small d at modest time 
scales. 
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(a) d=l 
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(b) d=b 
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(c) d=5 
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(d) d=10 

FIG. 8. Long-time behavior of the exact solution for Gaus¬ 
sian initial conditions with = 20, = 21,G OO - 1 • 



g G m? 

2320 1.3760 2.0000 0.9359 2.9359 

5628 1.0572 1.9805 2.1877 4.1739 

7266 1.0269 1.9315 3.1255 5.0638 

8728 1.0129 1.8596 4.4451 6.2765 

0000 1.0000 1.2500 ^/d^ 


TABLE I. Initial values for mr = 2.0, Qaa = 1. 


d 

9 

G 


r 

ro 

1 

2.4363 

1.3760 

21.0000 

0.9359 

2.9359 

5 

5.9798 

1.0572 

20.9805 

2.1877 

4.1865 

10 

7.9590 

1.0269 

20.9315 

3.1255 

5.1199 

20 

10.0209 

1.0129 

20.8596 

4.4450 

6.4305 

CX3 

20.0000 

1.0000 

20.2500 

Vd-1 



TABLE 11. Initial values for 






















